   This file documents the GNU Scientific Library.

   Copyright (C) 1996, 1997, 1998, 1999 The GSL Project.

   Permission is granted to make and distribute verbatim copies of this
manual provided the copyright notice and this permission notice are
preserved on all copies.

   Permission is granted to copy and distribute modified versions of
this manual under the conditions for verbatim copying, provided that
the entire resulting derived work is distributed under the terms of a
permission notice identical to this one.

   Permission is granted to copy and distribute translations of this
manual into another language, under the above conditions for modified
versions, except that this permission notice may be stated in a
translation approved by the Foundation.

GSL
***

   This file documents the GNU Scientific Library, a collection of
numerical routines for scientific computing.

   As of 19 May 2000 the library is in developers release only and is
not recommended for general use.

Preliminaries
*************

   The GNU Scientific Library (GSL) is a collection of routines for
numerical computing.  The routines are written from scratch by the GSL
team (see Contributors to GSL.) in C, and are meant to present a
modern Applications Programming Interface (API) for C programmers, while
allowing wrappers to be written for very high level languages.

   GSL is currently in developers release, for people who want to work
on the library itself. When the library is complete and fully tested it
will be announced for general use.

Using the library
*****************

   This chapter describes how to compile programs that use GSL, and
introduces its conventions.

ANSI C Compliance
=================

   The library is written in ANSI C and is intended to conform to the
ANSI C standard.  It should be portable to any system with a working
ANSI C compiler.

   The library does not rely on any non-ANSI extensions in the
interface it exports to the user.  Programs you write using GSL can be
ANSI compliant.  Extensions which can be used in a way compatible with
pure ANSI C are supported, however, via conditional compilation.  This
allows the library to take advantage of compiler extensions on those
platforms which support them.

   When an ANSI C feature is known to be broken on a particular system
the library will exclude any related functions at compile-time.  This
should make it impossible to link a program that would use these
functions and give incorrect results.

   To avoid namespace conflicts all exported function names and
variables have the prefix `gsl_', while exported macros have the prefix
`GSL_'.

Compiling and Linking
=====================

   The library header files are installed in their own `gsl' directory.
You should write any preprocessor include statements with a `gsl/'
directory prefix thus,

     #include <gsl/gsl_math.h>

If the directory is not installed on the standard search path of your
compiler you will also need to provide its location to the preprocessor
as a command line flag.  The default location of the `gsl' directory is
`/usr/local/include/gsl'.

   The library is installed as a single file, `libgsl.a'.  A shared
version of the library is also installed on systems that support shared
libraries.  The default location of these files is `/usr/local/lib'.
To link against the library you need to specify both the main library
and a supporting BLAS library.  A suitable blas implementation is
provided in `libgslblasnative' if your system does not provide one.
The following example shows how to link an application with the library,

     gcc app.o -lgsl -lgslblasnative -lm

The following command line shows how to link the same application with
an alternative blas library `libmycblas',

     gcc app.o -lgsl -lgslblascblas -lmycblas -lm

The library must conform to the CBLAS standard.  For more information
see BLAS Support.

   The program `gsl-config' provides information on the local version
of the library.  For example, the following command shows that the
library has been installed under the directory `/usr/local',

     bash$ gsl-config --prefix
     /usr/local

Further information is available using the command `gsl-config --help'.

Shared Libraries
================

   To run a program linked with the shared version of the library it
may be necessary to define the shell variable `LD_LIBRARY_PATH' to
include the directory where the library is installed.  For example,

     LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH ./app

To compile a statically linked version of the program instead, use the
`-static' flag in `gcc',

     gcc -static app.o -lgsl -lgslblasnative -lm

Automake macros
===============

   The GSL library also provides some useful automake macros to use in
your applications. The `gsl.m4' file contains the definitions needed to
detect GSL in automatic Makefile generation.

   To use the macros simply add the following line to your
`configure.in' file:

     AM_PATH_GSL(GSL_VERSION,AC_MSG_ERROR(Error message!))

where `GSL_VERSION' is the GSL Library version you need to compile and
run your application, and where "Error message" is the warning you want
to show if the required library is not found.

   Then, you can add the variables `GSL_LIBS' and `GSL_CFLAGS' to your
Makefile.am files to obtain the correct compiler flags.  `GSL_LIBS' is
equal to the output of the `gsl-config --libs' command and `GSL_CFLAGS'
is equal to `gsl-config --cflags' command. For example,

     libgsdv_la_LDFLAGS =    \
             $(GTK_LIBDIR) \
             $(GTK_LIBS) -lgsdvgsl $(GSL_LIBS) -lgslblasnative

Inline functions
================

The `inline' keyword is not part of ANSI C and the library does not
export any inline function definitions by default.  The inline versions
of functions can be included by defining the macro `HAVE_INLINE' when
compiling an application.

     gcc -c -DHAVE_INLINE app.c

If you use `autoconf' this macro can be defined automatically.  The
following test should be placed in your `configure.in' file,

     AC_C_INLINE

     if test "$ac_cv_c_inline" != no ; then
       AC_DEFINE(HAVE_INLINE,1)
       AC_SUBST(HAVE_INLINE)
     fi

and the macro will then be defined in the compilation flags or by
including the file `config.h' before any library headers.  If you do
not define the macro `HAVE_INLINE' then the slower non-inlined versions
of the functions will be used instead.

Long double
===========

   The extended numerical type `long double' is part of the ANSI C
standard and should be available in every modern compiler.  However, the
precision of `long double' is platform dependent, and this should be
considered when using it.  The IEEE standard only specifies the minimum
precision of extended precision numbers, while `double' is the same on
all platforms.

   In some system libraries the `stdio.h' formatted input/output
functions `printf' and `scanf' are not implemented correctly for `long
double'.  Undefined or incorrect results are avoided by testing these
functions during the `configure' stage of library compilation and
eliminating certain GSL functions which depend on them if necessary,

     checking whether printf/scanf works with long double... no

Consequently when `long double' formatted input/output does not work on
a given system it will not be possible to link a program which uses GSL
functions dependent on this.

   If it is necessary to work on a system which does not support
formatted `long double' input/output then the options are to use binary
formats or to convert `long double' results into `double' for reading
and writing.

Alternative optimized functions
===============================

   The main implementation of some functions in the library will not be
optimal on all architectures.  For example, there are several ways to
compute a Gaussian random variate and their relative speeds are
platform-dependent.  In cases like this the library provides alternate
implementations of these functions with the same interface.  If you
write your application using calls to the standard implementation you
can select an alternative version later via a preprocessor definition.
It is also possible to introduce your own optimized functions this way
while retaining portability.  For example,

     #ifdef SPARC
     #define gsl_ran_gaussian gsl_ran_gaussian_ratio_method
     #endif
     #ifdef INTEL
     #define gsl_ran_gaussian my_gaussian
     #endif

Note that the alternative implementations will not produce bit-for-bit
identical results, and in the case of random number distributions will
produce an entirely different stream of random variates.

Support for different numeric types
===================================

   Many functions in the library are defined for different numeric
types.  This feature is implemented by varying the name of the function
with a type-related modifier -- a primitive form of C++ templates.  The
modifier is inserted into the function name after the initial module
prefix.  The following table shows the function names defined for all
the numeric types of an imaginary module `gsl_foo' with function `fn',

     gsl_foo_fn               double
     gsl_foo_long_double_fn   long double
     gsl_foo_float_fn         float
     gsl_foo_long_fn          long
     gsl_foo_ulong_fn         unsigned long
     gsl_foo_int_fn           int
     gsl_foo_uint_fn          unsigned int
     gsl_foo_short_fn         short
     gsl_foo_ushort_fn        unsigned short
     gsl_foo_char_fn          char
     gsl_foo_uchar_fn         unsigned char

The normal numeric precision `double' is considered the default and
does not require a suffix.  For example, the function `gsl_stats_mean'
computes the mean of double precision numbers, while the function
`gsl_stats_int_mean' computes the mean of integers.

   A corresponding scheme is used for library defined types, such as
`gsl_vector' and `gsl_matrix'.  In this case the modifier is appended
to the type name.  For example, if a module defines a new
type-dependent struct or typedef `gsl_foo' it is modified for other
types in the following way,

     gsl_foo                  double
     gsl_foo_long_double      long double
     gsl_foo_float            float
     gsl_foo_long             long
     gsl_foo_ulong            unsigned long
     gsl_foo_int              int
     gsl_foo_uint             unsigned int
     gsl_foo_short            short
     gsl_foo_ushort           unsigned short
     gsl_foo_char             char
     gsl_foo_uchar            unsigned char

When a module contains type-dependent definitions the library provides
individual header files for each type.  The filenames are modified as
shown in the below.  For convenience the default header includes the
definitions for all the types.  To include only the double precision
header, or any other specific type, file use its individual filename.

     #include <gsl/gsl_foo.h>                  All types listed below
     #include <gsl/gsl_foo_double.h>           double
     #include <gsl/gsl_foo_long_double.h>      long double
     #include <gsl/gsl_foo_float.h>            float
     #include <gsl/gsl_foo_long.h>             long
     #include <gsl/gsl_foo_ulong.h>            unsigned long
     #include <gsl/gsl_foo_int.h>              int
     #include <gsl/gsl_foo_uint.h>             unsigned int
     #include <gsl/gsl_foo_short.h>            short
     #include <gsl/gsl_foo_ushort.h>           unsigned short
     #include <gsl/gsl_foo_char.h>             char
     #include <gsl/gsl_foo_uchar.h>            unsigned char

Compatibility with C++
======================

   The library header files automatically define functions to have
`extern "C"' linkage when included in C++ programs.

Aliasing of arrays
==================

   The library assumes that arrays, vectors and matrices passed as
arguments are not aliased and do not overlap. This allows the library to
use additional optimizations, and removes the need to handle overlapping
memory regions as a special case.  If overlapping arguments are used
then the results of such functions will be undefined.

Error handling in GSL
*********************

   This chapter describes the way that GSL functions report and handle
errors.  By examining the status information returned by every GSL
function you can determine whether it succeeded or failed, and if it
failed you can find out what the precise cause of failure was.  You can
also define your own error handling functions to modify the default
behavior of the library.

Error reporting
===============

   GSL follows the thread-safe error reporting conventions of the POSIX
Threads library.  Functions in GSL return a non-zero error code to
indicate an error and `0' to indicate success.

     int status = gsl_function(...)

     if (status) { /* an error occurred */
       .....       /* the value of status specifies the type of error */
     }

   GSL routines report an error whenever they cannot perform the task
requested of them.  For example, a root-finding function would return a
non-zero error code if could not converge to the requested accuracy, or
exceeded a limit on the number of iterations.  Situations like this are
a normal occurrence when using any mathematical library and you should
check the return status of the GSL functions that you call.

   Whenever a GSL routine reports an error the return value specifies
the type of error.  The return value is analogous to the value of the
variable `errno' in the C library.  However, the C library's `errno' is
a global variable, which is not thread-safe (There can be only one
instance of a global variable per program.  Different threads of
execution may overwrite `errno' simultaneously).  Returning the error
number directly avoids this problem.  The caller can examine the return
code and decide what action to take, including ignoring the error if it
is not considered serious.

   The error code numbers are defined in the file `gsl_errno.h'.  They
all have the prefix `GSL_' and expand to non-zero constant integer
values.  Many of the error codes use the same base name as a
corresponding error code in C library.  Here are some of the most common
error codes,

 - Macro: int GSL_EDOM
     Domain error; used by mathematical functions when an argument
     value does not fall into the domain over which the function is
     defined (like EDOM in the C library)

 - Macro: int GSL_ERANGE
     Range error; used by mathematical functions when the result value
     is not representable because of overflow or underflow (like ERANGE
     in the C library)

 - Macro: int GSL_ENOMEM
     No memory available.  The system cannot allocate more virtual
     memory because its capacity is full (like ENOMEM in the C
     library).  This error is reported when a GSL routine encounters
     problems when trying to allocate memory with `malloc'.

 - Macro: int GSL_EINVAL
     Invalid argument.  This is used to indicate various kinds of
     problems with passing the wrong argument to a library function
     (like EINVAL in the C library).

Here is an example of some code which checks the return value of a
function where an error might be reported,

     int status = gsl_fft_complex_radix2_forward (data, length);

     if (status) {
         if (status == GSL_EINVAL) {
            fprintf (stderr, "invalid argument, length=%d\n", length);
         } else {
            fprintf (stderr, "failed, gsl_errno=%d\n", status);
         }
         exit (-1);
     }

The function `gsl_fft_complex_radix2' only accepts integer lengths
which are a power of two.  If the variable `length' is not a power of
two then the call to the library function will return `GSL_EINVAL',
indicating that the length argument is invalid.  The `else' clause
catches any other possible errors.

Error handlers
==============

   In addition to reporting errors the library also provides an optional
error handler.  The error handler is called by library functions when
they are about to report an error (for example, just before they
return).  The purpose of the handler is to provide a function where a
breakpoint can be set when running under the debugger.

   The default behavior of the error handler is to print a short message
and call `abort()' whenever an error is reported by the library.  If a
library routine reports an error then the whole program will core-dump.
This is a safe default for lazy programmers who do not check the
return status of library routines (we don't encourage you to write
programs this way).  If you turn off the default error handler or
provide your own error handler then it is your responsibility to check
the return values of the GSL routines.

   All GSL error handlers have the type `gsl_error_handler_t', which is
defined in `gsl_errno.h',

 - Data Type: gsl_error_handler_t
     This is the type of GSL error handler functions.  An error handler
     will be passed three arguments, specifying the reason for the
     error, the source file in which it occurred, and the line number
     in that file.  The source file and line number are set at compile
     time using the `__FILE__' and `__LINE__' directives in the
     preprocessor.  An error handler function returns type `void'.
     Error handler functions should be defined like this,

          void HANDLER (const char * reason, const char * file, int line)

To request the use of your own error handler you need to call the
function `gsl_set_error_handler' which is also declared in
`gsl_errno.h',

 - Function: gsl_error_handler_t gsl_set_error_handler
          (gsl_error_handler_t NEW_HANDLER)
     This functions sets a new error handler, NEW_HANDLER, for the GSL
     library routines.  The previous handler is returned (so that you
     can restore it later).  Note that the pointer to a user defined
     error handler function is stored in a static variable, so there
     can only be one error handler per program.

          old_handler = gsl_set_error_handler (&my_error_handler);

          .....     /* code uses new handler */

          gsl_set_error_handler (old_handler) ; /* restore old handler */

     To use the default behavior (`abort' on error) set the error
     handler to `NULL',

          old_handler = gsl_set_error_handler (NULL);

Here is a skeleton outline of a program which defines its own error
handler.  Imagine that the program does interactive data analysis -
there is a main loop which reads commands from the user and calls
library routines with user-supplied arguments,

     #include <setjmp.h>
     #include <gsl/gsl_errno.h>

     jmp_buf main_loop;
     void my_error_handler (const char *reason, const char *file, int line);

     main ()
     {
        gsl_set_error_handler (&my_error_handler);

        while (1)
          {
            .... /* read command from user */

            if (setjmp (main_loop) == 0)
              {
                 .... /* call GSL routines requested by user */
              }
            else
              {
                 .... /* my_error_handler bailed out, GSL gave an error */
              }
          }
     }

     void
     my_error_handler (const char *reason, const char *file, int line)
     {
         fprintf (stderr, "GSL error: %s\n", reason);
         longjmp (main_loop);
     }

Before entering the interactive loop the program uses
`gsl_set_error_handler' to provide its own error handler
`my_error_handler' for GSL error reports.  After this point the
function `my_error_handler' will be invoked whenever an error is
reported by GSL. The new error handler prints the cause of the error
(the string `reason') and then does a non-local jump back to the main
loop.  This would allow the user to fix the command which caused the
error and try again.

Using GSL error reporting in your own functions
===============================================

   If you are writing numerical functions in a program which also uses
GSL code you may find it convenient to adopt the same error reporting
conventions as in the library.

   To report an error you need to call the function `gsl_error' with a
string describing the error and then return an appropriate error code
from `gsl_errno.h', or a special value, such as `NaN'.  For convenience
`gsl_errno.h' defines two macros to carry out these steps:

 - Macro: GSL_ERROR (REASON, GSL_ERRNO)
     This macro reports an error using the GSL conventions and returns a
     status value of `gsl_errno'.  It expands to the following code
     fragment,

          gsl_error (reason, __FILE__, __LINE__, gsl_errno) ;
          return gsl_errno ;

     The macro definition in `gsl_errno.h' actually wraps the code in a
     `do { ... } while (0)' block to prevent possible parsing problems.

   Here is an example of how the macro could be used to report that a
routine did not achieve a requested tolerance.  To report the error the
routine needs to return the error code `GSL_ETOL'.

     if (residual > tolerance)
       {
         GSL_ERROR("residual exceeds specified tolerance", GSL_ETOL) ;
       }

 - Macro: GSL_ERROR_RETURN (REASON, GSL_ERRNO, VALUE)
     This macro is the same as `GSL_ERROR' but returns a user-defined
     status value of VALUE instead of an error code.  It can be used for
     mathematical functions that return a floating point value.

   Here is an example where a function needs to return a `NaN' because
of a mathematical singularity,

     if (x == 0)
       {
         GSL_ERROR_RETURN("argument lies on singularity", GSL_ERANGE, NAN) ;
       }

Mathematical Functions
**********************

   This chapter describes basic mathematical functions.  Some of these
functions are present in system libraries, but the alternative versions
in the library can be used as a substitute when the system functions are
not available.

   The functions and macros are defined in the header file `gsl_math.h'.

Mathematical Constants
======================

   The library ensures that the standard BSD mathematical constants are
defined. For reference here is a list of the constants.

`M_E'
     The base of exponentials, e

`M_LOG2E'
     The base-2 logarithm of e, \log_2 (e)

`M_LOG10E'
     The base-10 logarithm of e, \log_10 (e)

`M_SQRT2'
     The square root of two, \sqrt 2

`M_SQRT1_2'
     The square root of one-half, \sqrt{1/2}

`M_SQRT3'
     The square root of three, \sqrt 3

`M_PI'
     The constant pi, \pi

`M_PI_2'
     pi divided by two, \pi/2

`M_PI_4'
     pi divided by four, \pi/4

`M_SQRTPI'
     The square root of pi, \sqrt\pi

`M_2_SQRTPI'
     Two divided by the square root of pi, 2/\sqrt\pi

`M_1_PI'
     The reciprocal of pi, 1/\pi

`M_2_PI'
     Twice the reciprocal of pi, 2/\pi

`M_LN10'
     The natural logarithm of ten, \ln(10)

`M_LN2'
     The natural logarith of two, \ln(2)

`M_LNPI'
     The natural logarithm of pi, \ln(\pi)

`M_EULER'
     Euler's constant, \gamma

Elementary Functions
====================

 - Function: double gsl_log1p (const double X)
     This function computes the value of \log(1+x) in a way that is
     accurate for small X. It provides an alternative to the system
     function `log1p(x)'.

 - Function: double gsl_hypot (const double X, const double Y)
     This function computes the value of \sqrt{x^2 + y^2} in a way that
     avoid overflow. It provides an alternative to the system function
     `hypot(x,y)'.

Testing the Sign of Numbers
===========================

 - Macro: GSL_SIGN (x)
     This macro returns the sign of X. It is defined as `((x) >= 0 ? 1
     : -1)'. Note that with this definition the sign of zero is positive
     (regardless of its IEEE sign bit).

Testing for Odd and Even Numbers
================================

 - Macro: GSL_IS_ODD (n)
     This macro evaluates to 1 if N is odd and 0 if N is even. The
     argument N must be of integer type.

 - Macro: GSL_IS_EVEN (n)
     This macro is the opposite of `GSL_IS_ODD(n)'. It evaluates to 1 if
     N is even and 0 if N is odd. The argument N must be of integer
     type.

Maximum and Minimum functions
=============================

 - Macro: GSL_MAX (a, b)
     This macro returns the maximum of A and B. It is defined as `((a)
     > (b) ? (a):(b))'.

 - Macro: GSL_MIN (a, b)
     This macro returns the minimum of A and B. It is defined as `((a)
     < (b) ? (a):(b))'.

 - Function: extern inline double GSL_MAX_DBL (double A, double B)
     This function returns the maximum of the double precision numbers
     A and B using an inline function. The use of a function allows for
     type checking of the arguments as an extra safety feature. On
     platforms where inline functions are not available the macro
     `GSL_MAX' will be automatically substituted.

 - Function: extern inline double GSL_MIN_DBL (double A, double B)
     This function returns the minimum of the double precision numbers
     A and B using an inline function. The use of a function allows for
     type checking of the arguments as an extra safety feature. On
     platforms where inline functions are not available the macro
     `GSL_MIN' will be automatically substituted.

 - Function: extern inline int GSL_MAX_INT (int A, int B)
 - Function: extern inline int GSL_MIN_INT (int A, int B)
     These functions return the maximum or minimum of the integers A
     and B using an inline function.  On platforms where inline
     functions are not available the macros `GSL_MAX' or `GSL_MIN' will
     be automatically substituted.

 - Function: extern inline long double GSL_MAX_LDBL (long double A,
          long double B)
 - Function: extern inline long double GSL_MIN_LDBL (long double A,
          long double B)
     These functions return the maximum or minimum of the long doubles A
     and B using an inline function.  On platforms where inline
     functions are not available the macros `GSL_MAX' or `GSL_MIN' will
     be automatically substituted.

Complex Numbers
***************

   The functions described in this chapter provide support for complex
numbers.  The algorithms take care to avoid unnecessary intermediate
underflows and overflows, allowing the functions to evaluated over the
as much of the complex plane as possible. (FIXME: this still needs to be
done for the csc,sec,cot,csch,sech,coth functions).

   For multiple-valued functions the branch cuts have been chosen to
follow the conventions of Abramowitz and Stegun in the `Handbook of
Mathematical Functions'. The functions return principal values which are
the same as those is GNU Calc, which in turn are the same as those in
`Common Lisp, The Language (Second Edition)' (1) and the HP-28/48
series of calculators.

   ---------- Footnotes ----------

   (1) The first edition of `Common Lisp, The Language' used different
definitions from the second edition.

Complex numbers
===============

   Complex numbers are represented using the type `gsl_complex'. The
internal representation of this type may vary across platforms and
should not be accessed directly. The functions and macros described
below allow complex numbers to be manipulated in a portable way

   For reference, the default form of the `gsl_complex' type is given
by the following struct,

     typedef struct
     {
       double dat[2];
     } gsl_complex ;

The real and imaginary part are stored in contiguous elements of a two
element array. This eliminates any padding between the real and
imaginary parts, `dat[0]' and `dat[1]', allowing the struct to be
mapped correctly onto packed complex arrays.

 - Function: gsl_complex gsl_complex_xy (double X, double Y)
     This function uses the cartesian components (X,Y) to return the
     complex number z = x + i y.

 - Function: gsl_complex gsl_complex_polar (double R, double THETA)
     This function returns the complex number z = r \exp(i \theta) = r
     (\cos(\theta) + i \sin(\theta)) from the polar representation
     (R,THETA).

 - Macro: GSL_REAL (Z)
 - Macro: GSL_IMAG (Z)
     These macros return the real and imaginary parts of the complex
     number Z.

 - Macro: GSL_SET_COMPLEX (ZP, X, Y)
     This macro uses the cartesian components (X,Y) to set the real and
     imaginary parts of the complex number pointed to by ZP.  For
     example,

          GSL_SET_COMPLEX(&z, 3, 4)

     sets Z to be 3 + 4i.

 - Macro: GSL_SET_REAL (ZP,X)
 - Macro: GSL_SET_IMAG (ZP,Y)
     These macros allow the real and imaginary parts of the complex
     number pointed to by ZP to be set independently.

Properties of complex numbers
=============================

 - Function: double gsl_complex_arg (gsl_complex Z)
     This function returns the argument of the complex number Z,
     \arg(z), where -\pi < \arg(z) <= \pi.

 - Function: double gsl_complex_abs (gsl_complex Z)
     This function returns the magnitude of the complex number Z, |z|.

 - Function: double gsl_complex_abs2 (gsl_complex Z)
     This function returns the squared magnitued of the complex number
     Z, |z|^2.

 - Function: double gsl_complex_logabs (gsl_complex Z)
     This function returns the natural logarithm of the magnitude of the
     complex number Z, \log|z|.  It allows an accurate evaluation of
     \log|z| when |z| is close to one. The direct evaluation of
     `log(gsl_complex_abs(z))' would lead to a loss of precision in
     this case.

Complex arithmetic operators
============================

 - Function: gsl_complex gsl_complex_add (gsl_complex A, gsl_complex B)
     This function returns the sum of the complex numbers A and B,
     z=a+b.

 - Function: gsl_complex gsl_complex_sub (gsl_complex A, gsl_complex B)
     This function returns the difference of the complex numbers A and
     B, z=a-b.

 - Function: gsl_complex gsl_complex_mul (gsl_complex A, gsl_complex B)
     This functions returns the product of the complex numbers A and B,
     z=ab.

 - Function: gsl_complex gsl_complex_div (gsl_complex A, gsl_complex B)
     This functions returns the quotient of the complex numbers A and
     B, z=a/b.

 - Function: gsl_complex gsl_complex_add_real (gsl_complex A, double X)
     This functions returns the sum of the complex number A and the
     real number X, z=a+x.

 - Function: gsl_complex gsl_complex_sub_real (gsl_complex A, double X)
     This functions returns the difference of the complex number A and
     the real number X, z=a-x.

 - Function: gsl_complex gsl_complex_mul_real (gsl_complex A, double B)
     This functions returns the product of the complex number A and the
     real number X, z=ax.

 - Function: gsl_complex gsl_complex_div_real (gsl_complex A, double B)
     This functions returns the quotient of the complex number A and the
     real number X, z=a/x.

 - Function: gsl_complex gsl_complex_add_imag (gsl_complex A, double Y)
     This function returns the sum of the complex number A and the
     imaginary number iY, z=a+iy.

 - Function: gsl_complex gsl_complex_sub_imag (gsl_complex A, double Y)
     This function returns the difference of the complex number A and
     the imaginary number iY, z=a-iy.

 - Function: gsl_complex gsl_complex_mul_imag (gsl_complex A, double Y)
     This function returns the product of the complex number A and the
     imaginary number iY, z=a*(iy).

 - Function: gsl_complex gsl_complex_div_imag (gsl_complex A, double Y)
     This function returns the quotient of the complex number A and the
     imaginary number iY, z=a/(iy).

 - Function: gsl_complex gsl_complex_conjugate (gsl_complex Z)
     This function returns the complex conjugate of the complex number
     Z, z^* = x - i y.

 - Function: gsl_complex gsl_complex_inverse (gsl_complex Z)
     This function returns the inverse, or reciprocal, of the complex
     number Z, 1/z = (x - i y)/(x^2 + y^2).

 - Function: gsl_complex gsl_complex_negative (gsl_complex Z)
     This function returns the negative of the complex number Z, -z =
     (-x) + i(-y).

Elementary Complex Functions
============================

 - Function: gsl_complex gsl_complex_sqrt (gsl_complex Z)
     This function returns the square root of the complex number Z,
     \sqrt z. The branch cut is the negative real axis. The result
     always lies in the right half of the complex plane.

 - Function: gsl_complex gsl_complex_sqrt_real (double x)
     This function returns the complex square root of the real number
     X, where X may be negative.

 - Function: gsl_complex gsl_complex_pow (gsl_complex Z, gsl_complex A)
     The function returns the complex number Z raised to the complex
     power A, z^a. This is computed as \exp(\log(z)*a) using complex
     logarithms and complex exponentials.

 - Function: gsl_complex gsl_complex_pow_real (gsl_complex Z, double A)
     This function returns the complex number Z raised to the real
     power B, z^b.

 - Function: gsl_complex gsl_complex_exp (gsl_complex Z)
     This function returns the complex exponential of the complex number
     Z, \exp(z).

 - Function: gsl_complex gsl_complex_log (gsl_complex Z)
     This function returns the complex natural logarithm (base e) of
     the complex number Z, \log(z).  The branch cut is the negative
     real axis.

 - Function: gsl_complex gsl_complex_log10 (gsl_complex Z)
     This function returns the complex base-10 logarithm of the complex
     number Z, \log_10 (z).

 - Function: gsl_complex gsl_complex_log_b (gsl_complex Z, gsl_complex
          B)
     This function returns the complex base-B logarithm of the complex
     number Z, \log_b(z). This quantity is computed as the ratio
     \log(z)/\log(b).

Complex Trigonometric Functions
===============================

 - Function: gsl_complex gsl_complex_sin (gsl_complex Z)
     This function returns the complex sine of the complex number Z,
     \sin(z) = (\exp(iz) - \exp(-iz))/(2i).

 - Function: gsl_complex gsl_complex_cos (gsl_complex Z)
     This function returns the complex cosine of the complex number Z,
     \cos(z) = (\exp(iz) + \exp(-iz))/2.

 - Function: gsl_complex gsl_complex_tan (gsl_complex Z)
     This function returns the complex tangent of the complex number Z,
     \tan(z) = \sin(z)/\cos(z).

 - Function: gsl_complex gsl_complex_sec (gsl_complex Z)
     This function returns the complex secant of the complex number Z,
     \sec(z) = 1/\cos(z).

 - Function: gsl_complex gsl_complex_csc (gsl_complex Z)
     This function returns the complex cosecant of the complex number Z,
     \csc(z) = 1/\sin(z).

 - Function: gsl_complex gsl_complex_cot (gsl_complex Z)
     This function returns the complex cotangent of the complex number
     Z, \cot(z) = 1/\tan(z).

Inverse Complex Trigonometric Functions
=======================================

 - Function: gsl_complex gsl_complex_arcsin (gsl_complex Z)
     This function returns the complex arcsine of the complex number Z,
     \arcsin(z). The branch cuts are on the real axis, less than -1 and
     greater than 1.

 - Function: gsl_complex gsl_complex_arcsin_real (double Z)
     This function returns the complex arcsine of the real number Z,
     \arcsin(z). For z between -1 and 1, the function returns a real
     value in the range (-\pi,\pi]. For z less than -1 the result has a
     real part of -\pi/2 and a positive imaginary part.  For z greater
     than 1 the result has a real part of \pi/2 and a negative
     imaginary part.

 - Function: gsl_complex gsl_complex_arccos (gsl_complex Z)
     This function returns the complex arccosine of the complex number
     Z, \arccos(z). The branch cuts are on the real axis, less than -1
     and greater than 1.

 - Function: gsl_complex gsl_complex_arccos_real (double Z)
     This function returns the complex arccosine of the real number Z,
     \arccos(z). For z between -1 and 1, the function returns a real
     value in the range [0,\pi]. For z less than -1 the result has a
     real part of \pi/2 and a negative imaginary part.  For z greater
     than 1 the result is purely imaginary and positive.

 - Function: gsl_complex gsl_complex_arctan (gsl_complex Z)
     This function returns the complex arctangent of the complex number
     Z, \arctan(z). The branch cuts are on the imaginary axis, below -i
     and above i.

 - Function: gsl_complex gsl_complex_arcsec (gsl_complex Z)
     This function returns the complex arcsecant of the complex number
     Z, \arcsec(z) = \arccos(1/z).

 - Function: gsl_complex gsl_complex_arcsec_real (double Z)
     This function returns the complex arcsecant of the real number Z,
     \arcsec(z) = \arccos(1/z).

 - Function: gsl_complex gsl_complex_arccsc (gsl_complex Z)
     This function returns the complex arccosecant of the complex
     number Z, \arccsc(z) = \arcsin(1/z).

 - Function: gsl_complex gsl_complex_arccsc_real (double Z)
     This function returns the complex arccosecant of the real number Z,
     \arccsc(z) = \arcsin(1/z).

 - Function: gsl_complex gsl_complex_arccot (gsl_complex Z)
     This function returns the complex arccotangent of the complex
     number Z, \arccot(z) = \arctan(1/z).

Complex Hyperbolic Functions
============================

 - Function: gsl_complex gsl_complex_sinh (gsl_complex Z)
     This function returns the complex hyperbolic sine of the complex
     number Z, \sinh(z) = (\exp(z) - \exp(-z))/2.

 - Function: gsl_complex gsl_complex_cosh (gsl_complex Z)
     This function returns the complex hyperbolic cosine of the complex
     number Z, \cosh(z) = (\exp(z) + \exp(-z))/2.

 - Function: gsl_complex gsl_complex_tanh (gsl_complex Z)
     This function returns the complex hyperbolic tangent of the
     complex number Z, \tanh(z) = \sinh(z)/\cosh(z).

 - Function: gsl_complex gsl_complex_sech (gsl_complex Z)
     This function returns the complex hyperbolic secant of the complex
     number Z, \sech(z) = 1/\cosh(z).

 - Function: gsl_complex gsl_complex_csch (gsl_complex Z)
     This function returns the complex hyperbolic cosecant of the
     complex number Z, \csch(z) = 1/\sinh(z).

 - Function: gsl_complex gsl_complex_coth (gsl_complex Z)
     This function returns the complex hyperbolic cotangent of the
     complex number Z, \coth(z) = 1/\tanh(z).

Inverse Complex Hyperbolic Functions
====================================

 - Function: gsl_complex gsl_complex_arcsinh (gsl_complex Z)
     This function returns the complex hyperbolic arcsine of the
     complex number Z, \arcsinh(z).  The branch cuts are on the
     imaginary axis, below -i and above i.

 - Function: gsl_complex gsl_complex_arccosh (gsl_complex Z)
     This function returns the complex hyperbolic arccosine of the
     complex number Z, \arccosh(z).  The branch cut is on the real less
     than 1.

 - Function: gsl_complex gsl_complex_arccosh_real (double Z)
     This function returns the complex hyperbolic arccosine of the real
     number Z, \arccosh(z).

 - Function: gsl_complex gsl_complex_arctanh (gsl_complex Z)
     This function returns the complex hyperbolic arctangent of the
     complex number Z, \arctanh(z).  The branch cuts are on the real
     axis, less than -1 and greater than 1.

 - Function: gsl_complex gsl_complex_arctanh_real (double Z)
     This function returns the complex hyperbolic arctangent of the real
     number Z, \arctanh(z).

 - Function: gsl_complex gsl_complex_arcsech (gsl_complex Z)
     This function returns the complex hyperbolic arcsecant of the
     complex number Z, \arcsech(z) = \arccosh(1/z).

 - Function: gsl_complex gsl_complex_arccsch (gsl_complex Z)
     This function returns the complex hyperbolic arccosecant of the
     complex number Z, \arccsch(z) = \arcsin(1/z).

 - Function: gsl_complex gsl_complex_arccoth (gsl_complex Z)
     This function returns the complex hyperbolic arccotangent of the
     complex number Z, \arccoth(z) = \arctanh(1/z).

References and Further Reading
==============================

   The implementations of the elementary and trigonometric functions are
based on the following papers,

     T. E. Hull and Thomas F. Fairgrieve and Ping Tak Peter Tang,
     "Implementing Complex Elementary Functions Using Exception
     Handling", `ACM Transactions on Mathematical Software', Volume 20
     (1994), pp 215-244, Corrigenda, p553

     T. E. Hull and Thomas F. Fairgrieve and Ping Tak Peter Tang,
     "Implementing the complex arcsin and arccosine functions using
     exception handling", `ACM Transactions on Mathematical Software',
     Volume 23 (1997) pp 299-335

The general formulas and details of branch cuts can be found in the
following books,

     Abramowitz and Stegun, `Handbook of Mathematical Functions',
     "Circular Functions in Terms of Real and Imaginary Parts", Formulas
     4.3.55-58, "Inverse Circular Functions in Terms of Real and
     Imaginary Parts", Formulas 4.4.37-39, "Hyperbolic Functions in
     Terms of Real and Imaginary Parts", Formulas 4.5.49-52, "Inverse
     Hyperbolic Functions - relation to Inverse Circular Functions",
     Formulas 4.6.14-19.

     Dave Gillespie, `Calc Manual', Free Software Foundation, ISBN
     1-882114-18-3

Vectors and Matrices
********************

   The functions described in this chapter provide a simple vector and
matrix interface to ordinary C arrays. The memory management of these
arrays is implemented using a single underlying type, known as a block.
By writing your functions in terms of vectors and matrices you can pass
a single structure containing both data and dimensions without needing
additional function arguments.  The structures are compatible with the
vector and matrix formats used by BLAS routines.

Data types
==========

   All the functions are available for each of the standard data-types.
The versions for `double' have the prefix `gsl_block', `gsl_vector' and
`gsl_matrix'.  Similarly the versions for single-precision `float'
arrays have the prefix `gsl_block_float', `gsl_vector_float' and
`gsl_matrix_float'.  The full list of available types is given below,

     gsl_block                       double
     gsl_block_float                 float
     gsl_block_long_double           long double
     gsl_block_int                   int
     gsl_block_uint                  unsigned int
     gsl_block_long                  long
     gsl_block_ulong                 unsigned long
     gsl_block_short                 short
     gsl_block_ushort                unsigned short
     gsl_block_char                  char
     gsl_block_uchar                 unsigned char
     gsl_block_complex               complex double
     gsl_block_complex_float         complex float
     gsl_block_complex_long_double   complex long double

Corresponding types exist for the `gsl_vector' and `gsl_matrix'
functions.